Astronomy & Astrophysics manuscript no. PlanetCO-Regaly-vlO 
July 20, 2010 



©ESO2010 



Detectability of giant planets in protoplanetary disks 

by CO emission lines 

Zs. Regaly', Zs. Sandor^, C. P. DuUemond^, and R. van Boekel^ 



o 

(N 

0^ 



6 



en 
> 

On 

(N 

1^ 

O 
O 



> 

X 



' Konkoly Observatory of the Hungarian Academy of Sciences, P.O. Box 67, H-1525 Budapest, Hungary 

e-mail: regaly@konkoly.hu 
- Junior Research Group, Max-Planck Institut fiir Astronomic, Konigstuhl 17, D-691 17 Heidelberg, Germany 
' Max-Planck-Institut fiir Astronomic, Konigstuhl 17, D-691 17 Heidelberg, Germany 



Received March 15, 2010; accepted July 6, 2010 



ABSTRACT 



Context. Planets are thought to form in protoplanetary accretion disks around young stars. Detecting a giant planet still embedded in 
a protoplanetary disk would be very important and give observational constraints on the planet-formation process. However, detecting 
these planets with the radial velocity technique is problematic owing to the strong stellar activity of these young objects. 
Aims. We intend to provide an indirect method to detect Jovian planets by studying near infrared emission spectra originating in 
the protoplanetary disks around T Tauri stars. Our idea is to investigate whether a massive planet could induce any observable eifect 
on the spectral lines emerging in the disks atmosphere. As a tracer molecule we propose CO, which is excited in the ro- vibrational 
fundamental band in the disk atmosphere to a distance of ~ 2-3 AU (depending on the stellar mass) where terrestrial planets are 
thought to form. 

Methods. We developed a semi-analytical model to calculate synthetic molecular spectral line profiles in a protoplanetary disk using 
a double layer disk model heated on the outside by irradiation by the central star and in the midplane by viscous dissipation due to 
accretion. 2D gas dynamics were incorporated in the calculation of synthetic spectral lines. The motions of gas parcels were calculated 
by the publicly available hydrodynamical code FARGO which was developed to study planet-disk interactions. 
Results. We demonstrate that a massive planet embedded in a protoplanetary disk strongly influences the originally circular Keplerian 
gas dynamics. The perturbed motion of the gas can be detected by comparing the CO line profiles in emission, which emerge from 
planet-bearing to those of planet-free disk models. The planet signal has two major characteristics: a permanent line profile asymmetry, 
and short timescale variability correlated with the orbital phase of the giant planet. We have found that the strength of the asymmetry 
depends on the physical parameters of the star-planet-disk system, such as the disk inclination angle, the planetary and stellar masses, 
the orbital distance, and the size of the disk inner cavity. The permanent line profile asymmetry is caused by a disk in an eccentric 
state in the gap opened by the giant planet. However, the variable component is a consequence of the local dynamical perturbation 
by the orbiting giant planet. We show that a forming giant planet, still embedded in the protoplanetary disk, can be detected using 
contemporary or future high-resolution near-IR spectrographs like VLT/CRIRES and ELT/METIS. 
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1. Introduction 

According to the general consensus, planets and planetary sys- 
tems form in circumstellar disks. These disks consist of gaseous 
and solid dust material. There are two major mechanisms pro- 
posed for the formation of giant planets. The first is the disk in- 
stability, which requires a gravitationally unstable disk in which 
giant planets form by direct collapse of the gas as a conse- 
quence of its self-gravity, originally suggested by Kuiper|p951| l; 



[Cameron ( |1978] l and more recently by [Boss] (2001 1. Later on 
the gas giant may collect dust, which after settling to its center, 
forms a solid core. The second mechanism is the core-accretion 
process ( [Bodenheimer & Pollack[[1986[ [Pollack et"aL][1996| l, 
which is the final stage of planet formation in the planetesi- 
mal hypothesis ( |Safronov| [T972 ). In this hypothesis dust coag- 
ulates first and forms planetesimals (meter to kilometer-sized 
objects). The subsequent growth of planetesimals consists of 
two stages, the runaw ay ([Wetherill & Stewart] 1989| l and the 
oligarchic growth (Kokubo & Ida 1998[ ), in which the formed 
bodies are massive enough to increase their masses by gravity- 
assisted collisional accretion processes. In this way planetary 



embryos and, through their consecutive collisions, protoplan- 
ets, terrestrial planets, and planetary cores of giant planets can 
be formed. During the gas accretion phase the precursor of the 
proto giant planet reaches a critical mass (where the planetary 
envelope contains more material than its core), and at about the 
order of lOM^ starts to contract, causing an increased accretion 
rate, which in turn raises radiative energy losses resulting in a 
runaway gas accretion. In this model the planetary core rapidly 
builds up a massive envelope of gas from its surrounding disk, 
in less than lOMyr, and eventually a gas giant forms. Because 
the disk dispersal time is presumably about 5 Myr (Haisch et al. 
[200 ![ [Hillenbrand[ [2005] l, this model should also incorporate 
planetary migration to explain giant planets observed well inside 
the snow line, see [Alibert et al. (2004) . 

Both planetary formation processes have their weak points. 
Disk instability assumes a very massive disk that is gravita- 
tionally unstable or only marginally stable. Another as yet un- 
solved problem is that during the gravitational collapse an ef- 
fective cooling mechanism should work to ensure the formation 
of gravitationally bound clumps. While the radiative cooling is 
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not effective enough and the recently invoked thermal convec- 
tion is also ineffective, the disk instability as a common planet- 
formation process is questionable ( (Klahr 2008 ). 

Regarding the giant planet formation in the planetesimal hy- 
pothesis, ^OTdasmTetaL (2009]l convincingly presented the suc- 
cess of the core accretion model to produce giant planets within 
the disk lifetime in a planet population synthesis model, but the 
artificial slow-down of the type I migration was necessary. Type I 
migration has such a short time scale ( lWard|1997| l that the planet 
inevitably will be engulfed by the host star within the disk life- 
time. Another unsolved issue of the planetesimal hypothesis is 
the so-called "meter-sized barrier problem": i.e. the formation of 
meter-sized bodies is impeded by their quick inward drift to the 
star ( Weidenschilling|1977 1, and by mutual disruptive collisions 
due to their high relative velocities above 1 m/s ( |Blum & Wurm| 
2008| l. If no other physical processes are taking place, these two 
mechanisms would impede the formation of the meter-sized and 
consequently the larger planetesimals. 

A wealth of information about the planet formation process 
would clearly be provided if newly formed planets, or at least 
their signatures, could be observed in still gas-rich protoplan- 
etary disks. There are already attempts to discover signatures 
of giant planet-formation assuming the disk instability mecha- 



nism in the works of Narayanan et al. ( 2006 1 and Jang-Condell 



& Boss (2007). In the first study, for instance, the authors sug 
gest observations of HCO^ emission lines as a tracer of the ac- 
cumulation of large clumps of cold material, arguing that with 
high-resolution millimeter and sub-millimeter interferometers 
the cold clumps can be directly imaged. The observability of 
an already formed planet embedded in a circumstellar disk has 
also been studied recently by , Wolf & D'Angelo, ( ,2005] l and Wolf 
|et al.| ( |2007| l. In these works the authors investigated whether 
the influence of a forming planet on the protoplanetary disk can 
be detected by analyzing the spectral energy distribution (SED) 
coming from the star-disk system. It was found that by studying 
only the SED it is impossible to infer the presence of an embed- 
ded planet. On the other hand, it was also concluded that the dust 
re-emission in the hot regions near the gap could certainly be de- 
tected and mapped by ALMA for the nearby (less than lOOpc in 
distance) and approximately face-on protoplanetary disks. 



Another interesting attempt was presented by Clarke & 
( 2003[l, where the authors have investigated the possi- 
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bility of detection of giant planets by excess CO overtone emis- 
sion originating from the planetary accretion inflow in FU Ori 
type objects. Nevertheless, the periodic line profile distortions 
are observable only in tight systems in which the giant planet 
is orbiting within 0.25 AU. Another possibility to detect em- 
bedded giant planets is based on the radial velocity measure- 
ments, which favor high-mass, short-period companions. The 
detectability limit of an embedded giant planet by radial velocity 
measurements is ~ 10 Mj or ~ 2Mj (Mj is the mass of Jupiter) 
orbiting at 1 AU or 0.5 AU, respectively, due to heavy variability 



in the optical spectra of young (< 10 Myr) stars (e.g., [Paulson & 
Yelda| ( |2006l l; |Prato et al.| ( |2008| l). Note that although there were 



attempts to detect planets by radial velocity measurements in 
young systems, no firm detection has yet been repeated see, e.g., 
Setiaw an et"aL]([2008l) for TWHyab, debated later by|Huelamo| 
etal.'(2008). The discovery of planetary systems by direct imag- 
ing nowadays is restricted to large separations > 10 - lOOAU 
( Veras et al.|2009| ) and favors higher mass planets. As of today 
nine exoplanetary systems have been directly imaged in the op- 
tical or near-infrared band, but only five of them are still em- 
bedded: GQ Lup ( [Neuhauser et al.|2005[ ) and CT Cha ( [Schmidt 



the companion masses are, because high as 21.5 Mj and 17 Mj, 
respectively; 2M1207 ( [Chauvin et al.|[2004[ |2005a | hosting a 
4Mj mass planet; UScoCTIO 108 ([ Kashyap et al.[|2008 ) host- 
ing a 14Mj mass planet; j6 Pic (Lagra nge et al.|2009a|b[ ) hosting 
an 8Mj mass plane t. AB Pic ( .Chauvin et al.[[2005b[), HR 8799 
( [Marois et al.|2008) (triple system ), SR 1845 ( [BiUer et al.|2'006) l 
and Fomalhaut (Kalas et al."2008) are mature planetary systems 
with ages of about 30 Myr, 60 Myr, 100 Myr, and 200 Myr. Note 
that centimeter wavelength radio observations of a very young 
(< 0.1 Myr) disk around HLTau, presented by Greaves et al. 
(2008 1, revealed the possibility of a 12 Jupiter mass giant planet 
being in formation at ~ 75 AU. 

We investigate a new possibility to detect young giant plan- 
ets embedded in protoplanetary disks of TTauri stars via line 
profile distortions in its near-IR spectra. Contrary to the studies 
mentioned above, instead of dealing with the disk's density per- 
turbations, we investigate the gas dynamics in the disk, particu- 
larly near the gap opened by the planet. We have found that there 
is a reasonable possibility to discover giant planets embedded in 
the protoplanetary disk of young stars with our approach. 



Our idea has been inspired by recent findings of Kley & 
Dirksen| (2006 ), who showed that a massive giant planet embed- 
ded in a protoplanetary disk produces a transition of the disk 
from the nearly circular state into an eccentric state. The mass 
limit for this transition is 3 Mj for a disk with v = 10"^ uni- 
form kinematic viscosity (measured in dimensionless units, see 
details in Sect. [3.1[ ). The most visible manifestation of the disk's 
eccentric state is that the outer rim of the gap opened by the giant 
planet has an elliptic shape. After performing a series of hydro- 
dynamical simulations with the code FARGO ( Masset 2000| l, we 
can also confirm the results of Kley & Dirksen '{ 2006i) even us- 



ing (T-type viscosity. Note that beside the elliptic geometry that 
is clearly visible in the disk surface-density distribution, the gi- 
ant planet also disturbs the motion of the gas parcels near the 
gap. Using CO as tracer molecule of gas, which is excited to 
~ 3 - 5 AU (depending on the stellar mass) in the disk atmo- 
sphere ( fNajita et al.|2007| l, we investigate how the CO emission 
line profiles at 4.7 /zm are distorted by the gas dynamics. With 
our semi-analytical synthetic spectral model we explore the in- 
fluence of the mass of the giant planet, the disk inclination angle, 
the mass of the hosting star, the orbital distance of planet, and fi- 
nally, the size of inner cavity on the line profile distortions. 

The paper is structured as follows: in Sect. 2 we review the 
basic physics of our synthetic spectral line models. In Sect. 3 
we present disk models and details of numerical simulations. 
Our results of the CO ro-vibrational line profile distortions, the 
detectability of an embedded giant planet, and its observability 
constraints are presented in Sect. 4. The paper closes with the 
discussion of the results and concluding remarks. 



2. Synthetic spectral line model 

To calculate the CO ro-vibrational spectra emerging from pro- 
toplanetary disks we developed a semi-analytical line spectral 
model. The thermodynamical model of the disk is based on the 
two-layer flared disk model, originally proposed by Chiang & 
IGoldreich (1997 1, which describes the heating by stellar irra- 



|et al.[[20()8] l, but they are rather hosting a stellar companion as 



diation in a disk with a flared geometry. In our model we use 
the flared disk approximation, and assume the disk to consist 
of two layers: an optically thick interior, producing continuum 
radiation and an optically thin atmosphere, producing line emis- 
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Fig. 1. Double-layer flared disk model used in our simulations, 
after |Chiang & Goldreich (1997i. In our model the following 
emission components are included: the optically thin atmosphere 
emission, with a temperature Tatm above the optically thick disk 
interior, with the temperature T-^^i; the continuum emission of 
the rounded-ofF disk inner rim, with a temperature Trim and the 
stellar continuum assumed to be a blackbody with a temperature 
r*. The accretion heating of disk interior is also accounted for. 



sion or absorption in the spectraFJThe boundary between the 
two layers lies where the dust becomes optically thick in the vi- 
sual wavelengths, i.e. at optical depth ry = 1 along the line of 
the incident stellar irradiation. The irradiation from the boundary 
layer formed at the inner edge of the disk ( [Popham et al.|1993] l 
is neglected though because its existence had not yet been con- 
firmed conclusively by observations yet. Nevertheless, the emis- 
sion from the inner disk rim is substantial, resulting in strong 
near-IR bump in Herbig Ae/Be ( Natta et al.||2001| l and weaker 
one in T Tauri SEDs ([M uzerolle et al.|2003| l. In addition it was 
found by Monnier et al. (2006) that the simple vertical wall as- 
sumption for the inner rim is inconsistent with the interfero- 
metric measurements. Indeed, the shape of the disk inner rim 
is presumably rounded-ofF, resulting in less inclination-angle- 
dependent near-IR excess emission (Isella & Natta"2005"). A 
schematic representation of the disk model is shown in FiglT] If 
the disk atmosphere is superheated by the stellar irradiation with 
respect to its interior, we may expect emission lines. Conversely, 
if the disk interior is heated above the temperature of the disk 
atmosphere by viscous dissipation, for example due to an abrupt 
increase in accretion rate, spectral lines are expected in absorp- 
tion. 

The vertical structure of a geometrically thin disk can be de- 
rived by considering vertical hydrostatic equilibrium ( Shakura] 
& Sunyaevl[T9 73 ). Ignoring any contribution from the gravita- 
tional force of the disk itself, the vertical density profile would 
be set by the equilibrium of gas pressure, the stellar gravita- 
tion and centrifugal force owing to orbital motion, resulting in 
a Gaussian vertical density profile. Because the thickness of the 
disk atmosphere in the double layer model is principally deter- 
mined by the grazing angle of the stellar irradiation, which is 
narrow {6{R) <; 1) in our computational domain {R < 5 AU), 
we assume that the disk atmosphere has uniform vertical density 
distribution. Note that the surface density in the disk atmosphere 



' We only consider the inner part of disk (R < 5 AU). A significant 
part of the fundamental band CO ro- vibrational emission arises in R = 
2-3 AU, however the line-to-continuum flux at 4.7 yum is sensitive to 
material that lies up to 5 AU in our models. In this way the outer part of 
disk, where the disk interior becomes optically thin to its own and even 
to the radiation from superheated atmosphere, is not considered. 



according to Eq. ( |D.2| i does not depend on the surface density 
(S(/?)) distribution by definition. 

The expected flux emerging from the double-layer disk is 
the result of the radiation of gas emission from the optically thin 
disk atmosphere above the continuum emission of dust in the op- 
tically thick disk interior and inner rim. Because the stellar con- 
tinuum contributes considerably to the line flux in the near-IR 
band compared to the disk continuum, especially in disks with 
an inner cavity, we have to also take into account the stellar con- 
tinuum in the total flux. In our model the stellar radiation is taken 
to be blackbody radiation of the stellar effective surface tempera- 
ture. Because the gas temperature is regulated by collisions with 
dust grains and stellar X-ray heating (which is significant for a 
young star), the gas and dust components may be thermally un- 
coupled in the tenuous di sk atmosphere below a crit ical density 



of«„ ~ 10"- 10"* cm'^ (Chiang & Goldreich 1997 1. As a resuh 



the gas temperature can be as high as ~ 
gion where the column density is below 



et al. 



2004 



4000 - 500 0K in a re- 
-- 10^' cm 2 (iGlassgold 

T 



I. Because the gas column density is ~ 10^^ cm 



in our model disk atmosphere, we adopt thermal coupling of 
dust and gas, i.e. T„(R) = T^iR) = T^t^^i„(R), where T.^i^^i„(R) 
is the atmospheric temperature defined by Eq. (|A.5jl. Although 
IGlassgold et al.| p004; and |Kamp & DuTlemond, ( ,2004 ) sug- 
gest the existence of an overheated layer above the superheated 
disk atmosphere, we did not consider this for simplicity's sake. 
In this model, assuming local thermodynamic equilibrium, the 
monochromatic intensity at frequency v emitted by gas parcels 
at a distance R to the star and an azimuthal angle (p observed at 
inclination angle / (/ = 0° meaning face on disk) can be given by 



I(v,R,4>,i) = B{v,TUR))e 



-t(v,R,M 



B(v,rat„(/?))(l-e-(^'«-^'''), 



(1) 



where t(v, R, (p, i) is the monochromatic optical depth of the disk 
atmosphere at a frequency v along the line of sight, B{v, Tint(R)) 
and B{v, T^t^(R)) are the Planck functions of the dust tempera- 
ture Tint{R) in the disk interior, and the gas temperature TatmC^) 
in disk atmosphere, respectively. Regarding the disk inner edge, 
which has a temperature ^^^(/Jo), it is assumed that it radiates 
as a blackbody {Inmiv) - B(v, TrimiRo))), see Appendix[B]for de- 
tails. The total flux emerging from the protoplanetary disk seen 
by inclination angle / at frequency v can thus be given as 

r*' r^" RdRdd) 

^di.sk(v,0 = J J /(v,/?,0,O-^cos(O + 



)Ro 
tnm(v,R ) 

D2 



D2 

Arini(0cOS(90-0, 



(2) 



where Rq and Ri are the disk inner and outer radii and D is the 
distance of the source to the observer. In Eq.l2]AriniCOs(90 - 
is the visible area of the disk inner rim. If the disk inner rim 
were a perfect vertical wall, its surface area would be A,in,(/) = 
Anh{R())R^y where /;(/?()) is the disk aspect ratio at the inner edge, 
which is taken to be 0.05. In order to be consistent with the ob- 
servations regarding the inclination-angle dependence of the rim 
flux ( |Isella & Natta|200 5^, the rim is taken to be a vertical wall 
seen under a constant 60 degree inclination angle. 

To calculate the synthetic spectra numerically, we first set 
up the two-dimensional computational domain with N^ loga- 
rithmically distributed radial and A^^ equidistantly distributed az- 
imuthal grid cells. For a given set of stellar parameters, shown in 
Table [T] first the grazing angle of the incident stellar irradiation 
is determined according to Eq. (A.2i. The temperature distribu- 
tions in the disk atmosphere and interior are calculated according 
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to Eq. ( |A3] i and Eqs. (|AJ]|X9||R2]|BJ]|, respectively. After de- 
termining the dust surface-density of the disk atmosphere with 
Eq. (D.2i, the optical depth along the line of sight is calculated 
by Eq. ( |D.4|i usi ng the gas opacity given by Eq. ( |E.1[ ) and consid- 
ering Eqs. (|E.4|)-(E.8 1 to calculate appropriate intrinsic line pro- 



files and Doppler shifts. In order to determine the mass absorp- 
tion coefficient of the '^C'^O molecule we used the transition 
data provided by Goorvitch (1994|l, such as transition probabil- 
ity A, lower and upper state energy E\, E^ and statistical weight 
^u- Knowing the optical depth (t(v, /?, 0, /)), the temperature in 
disk atmosphere and interior in the 2D computational domain 
the expected flux along the line of sight can be calculated by ap- 
plying Eq. ([T]i and Eq. Q at each frequency in the vicinity of the 
investigated transition. 

3. Disk models 

In this section we describe in detail our disk models that fed 
into the synthetic spectral line model. As a first step, we calcu- 
lated the disk surface-density and velocity distribution perturbed 
by a massive embedded planet with the publicly available hy- 
drodynamical code FARGO of Masset (2000). A very important 
difference to the planet-free case is that the disk surface-density 
distribution is very heavily modified. It shows a clear elliptic 
character of the gap opened by the giant planet. According to 
Kley & Dirksen| ( |2006[ ), |D'Angelo et al.|p006 ), and also to our 
calculations (see below), the disk becomes eccentric after sev- 
eral hundred orbits of the embedded planet. With regard to our 
synthetic spectral line model, the most relevant are those ef- 
fects that are originating from the gas dynamics caused by the 
giant planet. It is well known that a pure Keplerian motion of 
the gas parcels on circular orbits results in symmetric double- 
peaked emission line profiles ( [Home & Marsh] 1986| l. As the CO 
emission is strongly depressed 2-3 AU in our models, the ori- 
gin of the double-peaked profiles is clear Any deviation of the 
gas dynamics from the pure rotation will break this symmetry, 
resulting in the distortion of line profiles. If the distorting effect 
is strong enough, the line profile distortion could be significant 
to be detected by high-resolution spectroscopic observations. 

Indeed several TTauri stars show symmetric but centrally 
peaked CO line profiles that can be explained by radially extent 
CO emission overriding the double peaks ( [Najita et al.||2003| 
Brittain et al.||2009 ). Stronger CO emission can be expected if 
UV fluorescence plays a role (Krotkov et al.||1980|), or the gas 
temperature is not well coupled to the dust (Kamp & Dull emond| 



2004 Glassgold et al. 2004 1. In the tenuous region above the 



disk atmosphere the gas can be significantly hotter than is pre- 
dicted by the double-layer model. As a consequence the slowly 
rotating distant disk parcels could produce a substantial contri- 
bution to the low-velocity part of the line profile, resulting in a 
centrally peaked profile. Nevertheless, for simplicity's sake here 
we do not consider any of the above mentioned effects. 

3. 1 . Hydrodynamical setup 

We adopt dimensionless units in hydrodynamical simulations, 
for which the unit of length and mass is taken to be the orbital 
distance of the planet, and the mass of the central star, respec- 
tively. The unit of time, to, is taken to be the reciprocal of the 
orbital frequency of the planet, resulting in fo = 1 /27r, setting 
the gravitational constant to unity. In each disk model we as- 
signed four planetary masses to the embedded planet, which are 
in increasing order q = mpi/m. = 0.001, 0.003, 0.005, and 0.008, 
where nipi and m» are the planetary and stellar masses. Regarding 



the disks geometry, the flat thin disk approximation was assumed 
with an aspect ratio of H(R)/R - 0.05. The disk extends be- 
tween 0.2 - 5 dimensionless units. This computational domain 
was covered by 256 radial and 500 azimuthal grid cells. The ra- 
dial spacing was logarithmic, while the azimuthal spacing was 
equidistant. This results in practically quadratic grid cells, be- 
cause the approximation AR ~ RA(p is valid at each radius. The 
disks are driven by a-type viscosity ( Shakura & Sunyaev|19731 l, 
which is consistent with our thermal disk model with an inter- 
mediate value of a, which was taken to be 1 x 10"^. Note that 
assuming a cold disk (therefore virtually non-ionized) the vis- 
cosity generated by magnetohydrodynamical turbulence is hard 
to explain, although, according to Mukhopadhyay^ (200 8), the 
transient growth of two or three-dimensional pure hydrodynamic 
elliptic-type perturbations could result in a ^ 10"' - 10"^, de- 
pending on the disk scale height (decreasing a with increasing 
scale height). The disk's initial surface density profile is given 
by a power law 1.{R) - EqT?"'^^, where the surface density at 1 
distance unit Eq is taken to be 2.15 x 10"^ in dimensionless units. 
For simplicity, a locally isothermal equation of state is applied 
for the gas, and the disk self-gravity is neglected. While the inner 
boundary of the disk is taken to be open, allowing the disk mate- 
rial to leave the disk on accretion timescale, the outer boundary 
is closed, i.e. no mass supply is allowed. 

Below we assign physical units to the dimensionless quan- 
tities. If the stellar mass is assumed to be w, = IMq, the 
planetary masses in our models correspond to IMj, 3Mj, 5Mj, 
and 8Mj. Thanks to the dimensionless calculations it is possi- 
ble to scale the results to different stellar masses. For differ- 
ent stellar masses, the mass of the giant planet scales with the 
stellar mass, i.e mp\ = 1 x lO^^'m*, 3 x lO^^'m*, 5 x lO^^'m*, 
8 X 10"^ m», and the stellar masses are taken to be nit = 0.5 Mq, 
1 Mq, 1.5 M0. Furthermore, the distance unit is set to the dis- 
tance of the planet to the star, i.e., the planet orbits at 1 AU. 
Below we extended our models to tight and wide systems, where 
the giant planets are orbiting at 0.5 and 2 AU, respectively. In 
1 Mq stellar mass disk models the surface density at 1 AU is 
lo = 2.15 X IO^'^MqAU^^ ^ I9lg/cm2 at the beginning of 
simulation, which is also scaled with the mass of the central star. 
But note that the solution to the hydrodynamical system of equa- 
tions is independent of Eq, if there is no back-reaction to the 
planet, see the vertically integrated Navier-Stokes equations in 
■Kley(1999). For 0.5 Mq, 1 Mq and 1.5 Mq steUar mass model 
the disk mass in the computational domain (0.2 AU < R < 5 AU) 
is 5 X 10""* Mq, 1 X 10"^ Mq and 1.5 x 10"^ Mq. These values 
refer to the mass of the inner disk only, as the whole disk may 
extend to several tens or hundreds of AU, and the total disk mass 
is in the conventional range of 0.01 - 0. 1 Mq. More details about 
the models can be found in Table [U 



3.2. Synthetic spectral calculation setup 

We calculate the spectral line of a fundamental band CO transi- 
tion (V = 1 ^ OP(IO), /!() = 4.7545 jum), which is not blended 
by the higher excitation V=2— »1 and V=3— >2 lines, although 
in our thermal model the temperature is not high enough to ex- 
cite the higher vibrational levels at all |^ The CO rotational and 
vibrational level populations were calculated in local thermody- 
namical equilibrium. The line profiles were calculated in 200 
wavelengths and in the same numerical domain that was used in 



' In our thermal model the temperature stayed below 1500K every- 
where in the computational domain. At this temperature the V > 2 vi- 
brational levels are not thermally excited. 
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Table 1. Stellar and disk parameters in our hydrodynamical 
models. The stellar parameters were taken from a publicly avail- 
able tool ( |Siessetal.|2000[ ). The disk extends between 0.2-5 AU 
and the orbital distance of the massive planets is 1 AU. 



Model No. 


m, (Mq) 


r.(K) 


R. (Re) 


in 


mp (Mj) 


#1 


0.5 


3760 


1.4 


20,40,60 


0.5 


#2 


0.5 


3760 


1.4 


20,40,60 


1.5 


#3 


0.5 


3760 


1.4 


20,40,60 


2.5 


#4 


0.5 


3760 


1.4 


20,40,60 


4 


#5 


1 


4266 


1.83 


20,40,60 


1 


#6 


1 


4266 


1.83 


20,40,60 


3 


#7 


1 


4266 


1.83 


20,40,60 


5 


#8 


1 


4266 


1.83 


20,40,60 


8 


#9 


1.5 


4584 


2.22 


20,40,60 


1.5 


#10 


1.5 


4584 


2.22 


20,40,60 


4.5 


#11 


1.5 


4584 


2.22 


20,40,60 


7.5 
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hydrodynamical simulations. According to our preliminary tests 
on "circularly Keplerian'rldisks, the CO fundamental band spec- 
tra did not show significant changes with increasing outer disk 
boundary. This can be explained by the fact that a substantial 
part of the CO fundamental band emission is arising in the in- 
ner disk, consequently the outer part of the disk {R > 3 - 5 AU, 
depending on the stellar mass) does not contribute to the disk 
emission at 4.7 micron. Contrary to this, the distance of the inner 
boundary of the disk to the central star heavily influences the CO 
fundamental band spectra, i.e. the closer the disk inner boundary 
to the star the line-over-continuum is weaker and broader. This 
is because the disk innermost regions give stronger contribution 
to the continuum, which in turn weakens the CO emission nor- 
mahzed to continuum, although the CO emission itself is also 
getting stronger due to an increased amount of hot CO. The line 
broadening is a natural consequence of the larger orbital velocity 
of gas orbiting closer to the star 

Because we investigate disks with akeady formed planets, it 
is reasonable to assume that a cavity is formed at the very inner 
disk like the ones found in disks of several T Tauri stars ( Akeson| 
erar.^2005. Eisner et al.|2005[ [20071 |Salyk et al._2009 ). If dust 



evaporation is takes place close to the star, the CO will be de- 
pleted too, due to the disappearance of the dust which protects 
the CO molecules against UV photo dissociation. In this way, the 
observed inner radius of CO (^co) should intuitively be similar 
or slightly larger than the dust sublimation radius R^ub- Indeed, 
Rco was found to be larger than R^uh by a factor of a few in sev- 
eral disks, e.g., see Fig. 5. of Eisner et al. (2007 1 or Fig. 10 in 
|Salyketal.| ( [2509l ). Contrary to this, |Carr| 
CO inner radius is ~ 0.7 that of the dust. If the gas disk extends 
as far as 0.05 AU, there would be an additional broad component 
to the emission that is not modeled here. Taking these argumen- 
tations into account it is reasonable to assume that the CO inner 
radius is set by the dust sublimation. The dust sublimation radius 
/?sub is given by 
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Fig. 2. Logarithmic density {green shaded in gjcnr') and tem- 
perature contours (red dashed contours in K) of the disk. The 
density contours are stopped at a gas density p - 10""''* g/cm-' to 
avoid color crowding, while the maximum is p = 10"" g/cnr'. 
The height of the disk atmosphere, i.e the range where the ra- 
dial optical depth is less than 1 is shown with white dotted lines 
in the disk inner region. Comparing the unperturbed case (top) 
to the disk in which a gap exists at 0.5 AU < R < 2 AU, after 
depleting the density by 1/1000 (bottom), it is appreciable that 
however the disk atmosphere is decreased in height, the temper- 
ature distribution is not considerably changed. 



( 2007| l found that the where Eq. (|A.5|l with standard dust opacity law, /3 = 1 ([Rodmann 



R. 



0.03 



r. 



5/2 



R* 



AU, 



(3) 



We note that the term "Keplerian", widely used in the literature, 
here means "circularly Keplerian" motion. On the other hand, elliptic 
motion is also Keplerian, but below we will use the widely accepted 
nomenclature for circular motion, though we think it is not entirely cor- 
rect. 



et al. 2006 1 was used. Assuming T^ub - 1500 K for dust subUma- 
tion temperature, r» ~ 4000 K and Rt = I.^Rq for stellar sur- 
face effective temperature and radius, appropriate for a 2.5 Myr 
solar mass star, we obtain /?sub - 0.05 AU. To keep this simple, 
we set the disk inner boundary (namely the CO inner radius) 
fixed to ^co = 4/?sub = 0.2 AU, which is an acceptable value in 
disks hosted by T Tauri stars with 3500-5000 K surface tempera- 
ture and 1 .4 - 2.2 Rq radius. Note that /?sub was fixed throughout 
our models (except the ones where we investigated the effect of 
the size of the inner cavity) to make them comparable, neglect- 
ing the efifect of the change in stellar luminosity on the R^^b- 

It is reasonable to assume that a disk that is 2.5 Myr old con- 
tains a considerable amount of gas and tracer CO also in the 
inner disk, hence the stellar age is taken to be 2.5 Myr in all 
models. Note that according to |Haisch et al.| ( [2001| l, half of the 
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observed young stars in nearby embedded clusters (NGC 2024, 
Trapezium, IC 348, NGC 2264, NGC 2362, NGC 1960), with 
ages of about 3 Myr, show near-IR excess, which is related to 
the excess emission of dust in disks. Furthermore, Hillenbrand] 
( 2005 | l also concluded that the median lifetime of the optically 
thick inner disk is between 2-3 Myr. 

In order to simplify the models it is assumed that the dust 



consists of pure silicates with 0.1 /im grain size (Draine & Lee 



1984 1. According to this the mass absorption coefficient of the 
dust is taken to be 2320cm~/g at visual wavelength. Note that 
neither the effect of the coagulation nor the settling out of tenu- 
ous atmosphere of grains with a large size are taken into account. 
The size distribution and therefore the overall mass absorption 
coefficient of the dust is taken to be the same in the interior and 
the atmosphere. In this way the dust opacity is slightly overesti- 
mated. According to Eq. (D.4i the atmospheric gas density and 
thus gas line emission strength is slightly underestimated. The 
dust-to-gas and CO-to-gas mass ratios are assumed to be con- 
stant throughout the disk and are taken to be 10"^ and 4 x 10 "* 
(measured per gram of gas+dust mass) in the disk atmosphere 
and interior, respectively. 

For the planet-free disks circular Keplerian velocity distribu- 
tion is applied to determine the observed line center shift, given 
by Eq. (E.8i. For a massive planet-bearing disk observed un- 



der inclination angle /, the Doppler shift of the emission from 
a single patch of gas in the disk compared to its fundamental 
frequency vo emerging from a given R, point is calculated by 



Av(R, (j), i) 



V{) 



u^(R, (p) [sin(<^) + cos(0)] -H 
^(R, (j)) [cos(0) - sin(0)]} sin(0, 



(4) 



where the radial u^(R, (p) and azimuthal u^iR, (p) velocity com- 
ponents of gas parcels are provided by the hydrodynamic simu- 
lations. 

The perturbations in the surface density distribution of the 
disk interior and atmosphere were not taken into account, be- 
cause the disk atmosphere density given by Eq. (D.2i is inde- 
pendent of the density distribution of the disk even in the gap, 
as long as there is enough dust in the gap to keep it optically 
thick. To test the optical thick assumptions, let us for a mo- 
ment neglect the effect of disk shelf-shadowing that occurres 
in the inner edge of the gap. Applying the dust density in the 
disk I.i(R = 1 AU) - I.()Xd, where X^ = 0.01 is the dust-to- 
gas ratio and the dust opacity Ky - 2320cm^/g, in the gap, 
depleted to 0.1% of the surrounding density, the optical depth 
at visual wavelengths along the incident stellar irradiation is 
TV = I.d{R = 1 AlJ)Ky/6{R = 1 AU) ^ 1000. To argue that 
the shelf shadowing can be neglected we calculated the height of 
the disk atmosphere and temperature distribution in a planet-free 
and planet-bearing disk by a 2-D Monte-Carlo code RADMC 
for dust continuum radiative transfer ( |Dullemond & Dominik| 
|2004 ) with the same dust prescription as in our synthetic spec- 
tral model. In the planet-bearing disk model the gap is taken to be 
extended from 0.5 AU to 2 AU and artificially depleted in density 
to 0. 1 % that of the planet-free case. The amount of gap depletion 
is an average value measured in our hydrodynamic calculations. 
The results are shown in Fig. |2] It is evident that although the 
height of the disk atmosphere (ry - 1 along the line of sight 
of stellar irradiation) above the mid-plane is decreased in the 
gap, the temperature in the disk atmosphere is not substantially 
changed. Because the temperature is not significantly decreased 
we neglect the surface density perturbations. Note that the ef- 
fect of the distance of the gap from the central star on the shelf 



shadowing is not considered. The dust and gas temperature de- 
coupling was also not taken into account. Because the gas has a 
temperature well in excess of the d ust in depleted regions , where 
the gas column density is « 10^^ ( Glassgold et al.|2004[ Kamp 
|& DuUemond 20041), we can expect increased contribution to 
the CO line flux originating from the gap. Contrary to this, if the 
base of the gap is shadowed completely by the inner wall, the gap 
does not contribute to the CO emission. Thus a gap could cause 
a strong permanent line profile asymmetry due to its asymmetric 
geometry (see Fig. [5]), which requires further investigation. 

Because the accretion heating has been taken into account 
according to Eq. ( |A.9| l, we had to set an appropriate accretion 
rate. The accretion rate measured in hydrodynamical simula- 
tions gives a value of about 2 x 10"^ Mq/yi. Note that according 
to |Fang et al. ( 2009| l the accretion rate inferred from Ha emis- 
sion luminosity in the young star population of Lynds 1630N and 
1641 clouds in the Orion GMC with an age about 2-3 Myr, is 
between 10"'" - 10"*^ MQ/yr. In order to be consistent with the 
latter, the accretion rate is taken to be 5 x 10"'' Mq/yt for all our 
models. However, note that at this accretion rate the contribution 
of the viscous dissipation to the disk continuum is weak in the 
near-IR band. 



4. Results 

4.1. Effect of massive planets on the disk structure 

Three groups of models were computed through 2000 planetary 
orbits, with stellar and planetary parameters listed in Table [T] In 
each run the planet was kept fixed on circular orbit during the 
first 1000 orbits. After the 1000th orbit the planet was released. 
It thus felt the (gravitational) backreaction of the disk, resulting 
in its inward migration. In a good accordance with the expec- 
tations, we found that a more massive planet opened a broader 
and deeper gap. The depletion is 0.25%-0. 1% for planets with 
mass in the range of 1 Mj - 8 Mj. After the first 1000 orbits a 
quasi steady state disk structure has developed in all simulations, 
which slightly changed during the following 1000 orbits, when 
the planet orbit was not fixed anymore. 

To shed light on how the giant planet distorts the originally 
circularly Keplerian gas flow, several snapshots of density and 
velocity distributions were taken during the 2000th planetary or- 
bit. Figure[3] shows four snapshots of the density and the radial 
velocity component of the orbital velocity of gas for an 8 Mj 
mass planet orbiting an 1 Mq star (model #8) during one orbit. 
It is evident that in a disk with an embedded massive planet the 
overall orbits of gas parcels are noncircularly Keplerian because 
the radial components of their orbital velocity distribution are 
strongly departing from zero. Moreover, we found that the ellip- 
tic gap, while preserving its shape, precesses slowly retrograde 
with a period of about 150 planetary orbits. 

As we already mentioned, Kle y & Dirksen| ( [2006[ ) found that 
an originally circular disk with an embedded giant planet can 
reach an eccentric equilibrium state if the mass of the embed- 
ded planet is larger than a certain limit (3 Mj for a 1 Mq star). 
Kley & Dirksen found that for sufficiently wide gaps the growth 
of the eccentricity is induced by the interaction of the planet's 
gravitational potential with the disks material at the radial loca- 
tion of the 1 :3 (outer) Lindblad resonance. For smaller planetary 
masses, this effect is damped mainly by the co-orbital and the 1 :2 
Lindblad resonances. If the gap is deep and wide enough, which 
is the case for a giant planet, the above resonances cannot damp 
the eccentricity-exciting effect appearing at the 1:3 Lindblad res- 
onance, and the disk becomes eccentric. For a more detailed ex- 
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Fig. 3. Surface density {left four figures, orange colors in the range of 0.2 - 250g/cm^) and the radial velocity component of the 
orbital velocity in the disk plane {right four figures, blue colors in the range of -10.86 — i-10.96km/s) distributions of the disks 
gaseous material taken from four different azimuthal positions of the planet's 2000th orbit. In our model the total disk mass is 
0.05 Mq, the planet mass is 8 Mj, while the stellar mass is 1 Mq. The planetary orbit is shown with a white and red dashed circle 
in the density and velocity plots. The gap {in dark) and the elliptic shape of its outer rim is clearly visible in the surface density 
distribution. In the vicinity of the gap the disk material shows a strong deviation from the circular Keplerian rotation, see the white 
clumps in the radial velocity component distribution. 



planation of disk eccentricity growth see'Lubow(199r) and Kley] 
|& Dirksen, (.2006). The eccentric state of the disk in our simula- 
tions can also be clearly seen in Fig. [3] where the shape of the 
outer rim of the gap becomes elliptic in the surface density dis- 
tribution. 

The departure of the velocities of the gas parcels from the 
pure circularly Keplerian circular revolution can be character- 
ized by calculating their eccentricities. Considering the eccen- 
tric equilibrium state of the disk, we expect that each gas parcel 
would move on a non-circular orbit, which may be characterized 
most conveniently by an average eccentricity value. Therefore, 
for each disk radius between 0.2 AU and 5 AU, the azimuthally 
averaged eccentricities of the orbits of the gas parcels were cal- 
culated in the following way 



e{R) = I ^l + 2h{R, (p)c{R, cpfdcp. 

In Eq. (|5]l c{R, cj)) and h{R, tp) stand for 
c{R, (P) = x{R, cf>)vy{R, 0) - y{R, (/>)v,{R, 0) 
and 

v,{R,4>f + Vy{R,^f 



h{R, (f>) = 



1 



(5) 



(6) 



(7) 




2 y/x{R,(P)^+y{R,(l>)^' 

where v^{R,<p), Vy{R,<p) and x{R,(p), y{R,(p) are the Cartesian 
velocity components and coordinates at point /?, 0. We found 
that the averaged eccentricities differ considerably from zero, 
and each eccentricity curve reaches its maximum near the outer 
boundary of the gap (Fig.pl. A priori one would expect that the 
averaged eccentricity values are higher for more massive planets. 



R (AU) 

Fig. 4. Azimuthally averaged eccentricities as the functions of 
radii after 2000 orbits of the giant planet. The eccentricity curves 
in models #5-#8 (where the planetary masses are 1 Mj, 3 Mj, 
5 Mj and 8 Mj for 1 Mq stellar mass) are shown with dot-dashed, 
dashed, dotted, and solid lines, respectively. 



Indeed we found that the maximum of the azimuthally averaged 
eccentricity is monotonically increasing with a planetary mass 
in the range of 1 Mj < nipi < 8 Mj. The eccentricity curve for the 
8 Mj mass planet peaks about Cmax ~ 0.3, while the peak stays 
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Fig. 5. Offset of the line of sight component of the velocity dis- 
tribution from circularly Keplerian values on the disk surface in 
the hydrodynamic simulations, in the same model presented in 
Fig. [3] Snapshots are calculated for four different azimuthal po- 
sitions of the planet of the 2000th orbit. The velocity difference 
is in the range of -7.92 - 4.88 km/s. Here the inclination angle 
is taken to be / = 40° and the viewing angle is taken to be -90°, 
i.e the disk is seen from the bottom of Fig. [3] The planetary orbit 
is shown with a red dashed circle. Two permanent high-velocity 
regions can be seen at position angle 0° and 180°, and a variable 
pattern in the vicinity of the planet orbits. The strength of the 
variable component can reach that of the permanent one, e.g, in 
the panel at the bottom left (planet is at position angle 180°). 



well bellow 0. 1 for a 1 Mj mass planet. Note that a very similar 
behavior of the disk eccentricity has already been found by Kley 
|& Dirksen| ( |2006| ). In their cases, however, the maximum of the 
eccentricity curves is somewhat lower than in our cases, and at 
least 3 Mj is required to set the disk into eccentric state. This can 
be the consequence because contrary to Kley & Dirksen ( |2006 1, 
we used an a-type viscosity in our simulations. The higher ec- 
centricities we found are plausible inasmuch as |Kley & Dirksen 
( |200 6) found that the eccentricity of the disk is increasing with 
decreasing viscosity and in our approach the kinematic viscosity 
v{R) - aH^Q.]<^(R) measured in dimensionless units is 2. 5 x 10"^ 
atR - I, which is smaller than the v = 1 x 10"^ used by Kley & 
|Dirksen| ( |2006l ). 

As one can see, a giant planet has substantial impact on the 
density and velocity distributions of its host disk. It is an essen- 
tial question, whether the velocity perturbations appear in the 
line-of-sight velocity with substantial strength. FigurelS] shows 
the subtracted distribution of the line-of-sight velocity in the per- 
turbed and circularly Keplerian case. It is evident that the line- 
of-sight velocities show a significant departure from the circu- 
larly Keplerian fashion because the difference is non-vanishing. 
The radial velocity component of the orbital velocity and more 
importantly the line-of-sight velocity distributions have variable 
patterns following the planet on the top of a permanent excess 
seen at 0° and 180° position angle. Note that the the permanent 




Fig. 6. Asymmetric CO ro-vibrational line profile in a circularly 
Keplerian disk with a gap between 0.8 -2.8 AU (similar to model 
#8), where the gas is flowing in eccentric orbit with e = 0.2. The 
gas outside the gap is in Keplerian circular orbit. It is evident that 
the line profiles are asymmetric. For comparison, we display the 
symmetric double-peaked line profiles (dashed lines) emerging 
from a planet-free disk as well. 



pattern precesses slowly (with ~ 150 orbital periods), retrograde 
to the planet. The disk inclination angle / is taken to be 40° in 
the calculation of Fig. [5] and the disk was rotated to the line- 
of-sight in a way that the line-of-sight velocity component of 
dynamically perturbed gas parcels is maximized. Consequently 
we had expected that not only significant distortions appear in 
the line profiles, but that they vary in time within the orbital time 
scale of the giant planet. Because the deviation from circularly 
Keplerian velocity is in the range of -7.92-4.88 km/s, the width 
of variable component in the line profile should be ~ lOkm/s, 
depending on the inclination angle. 

4.2. Distortion of CO lines 

Below we address the following questions: (i) Does the line-of- 
sight component of the non-circularly Keplerian velocity distri- 
bution (Fig.Bll result in significant distortions in the CO spectral 
line profiles? (ii) Are these distortions varying in the planets or- 
bital timescale? (iii) How do the distortions depend on the incli- 
nation angle, the planetary and stellar mass, the orbital distance 
of the planet, the size of the inner cavity and the disk geome- 
try (if the disk is flared or not)? (iv) Are these distortions strong 
enough to be detected by high-resolution near-IR spectral mea- 
surements? A complete set of answers to the above questions 
may give us a novel method to detect massive planets embedded 
in protoplanetary accretion disks. 

As was shown in Fig.fflthe disk eccentricity is considerable 
near the gap in all planet-bearing disk models. If the gas parcels 
move on pure elliptic orbits in the gap, the CO line profile be- 
comes permanently asymmetric. This effect is illustrated in Fig. 
|6] where we have calculated the emerging line profile of a disk 
with an eccentricity e = 0.2 in the gap between 0.8 AU - 2.8 AU 
hosted by a 1 Mq star viewed 20°, 40° and 60° degree of inclina- 
tion angles. Thus we can expect that the non-circularly Keplerian 
property of the velocity field shown in Figs. l3] and l5] will break 
the symmetry of the spectral lines. Although the orbits of the gas 
parcels can be characterized by average eccentricities, their or- 
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bits are more complicated than those of a pure elhptic motion, 
because of the perturbations induced by the planet. 

First we examined the CO emission line profile distortions 
from a disk in which a 8Mj planet orbits a IMq star (model 
#8). As expected, the strongly non-circularly Keplerian veloc- 
ity flow significantly modified the CO line profiles. In Fig.lSta) 
we display the V=1-0P(10) line profiles for three different in- 
clination angles. The line profiles show a strongly asymmetric 
double-peaked shape. Studying Fig.l8](a) two distortion patterns 
of the spectral line profiles can be recognized. An excess can be 
seen near the red peak of the lines at about 11 km/s, 21 km/s 
and 29 km/s, for inclinations / = 20°, 40° and 60°, respectively. 
Moreover, a deficiency appears around the blue peak of each 
spectral line at slightly lower velocities. The above velocities 
correspond to the orbital velocity of the gas revolving around a 
1 Mq star at 0.8-0.9 AU, i.e. near the inner boundary of the gap, 
taking into account the appropriate inclination angles. Smaller- 
scale distortions also appear in the red and blue peaks as well, 
corresponding to the 0.8-2.9 AU region, i.e. the whole gap. 

Calculating the line profiles at different orbital phases dur- 
ing one orbit, we found that the shape of the already distorted 
line profile varies in time, yielding a clear dependence on the 
orbital position of the planet, see results in Fig. l9]for model #8. 
The permanent asymmetry in lines is caused by the permanent 
velocity pattern seen at PA 0° and 180° in Fig.|5] Note that the 
expression "permanent" is not entirely correct because the ve- 
locity pattern slowly (~ 150 orbital periods of planet) precesses 
retrograde to the planet, but we still use it henceforward. The 
variable component indicated with arrows is moving in regions 
between ~ ±25 km/s due to the variable pattern following the 
planet. Note that the maximum width of variable component 
is ~ 10 km/s, as is expected, which is resolvable by CRIRES, 
whose maximum resolution is 3 km/s. In this particular model 
the time scale of variations is on the order of weeks, because 
the orbital period of the planet is one year, thus the time elapsed 
between snapshots is approximately 18 days. 

4.3. Planetary and stellar masses 

To gain further insights, it is useful to study the influence of the 
masses on the strength of the spectral line distortions. First we 
have investigated the cases in which a 1 Mq star hosts a 1 Mj, 
3 Mj, 5 Mj, and 8 Mj mass planet at 1 AU, corresponding to mod- 
els #5-#8. Parts of the results are presented in Fig.lHJb) for model 
#5 (mp, - 1 Mj, m, = 1 Mq) and in Fig.lHta) and Fig.l9]for model 
#8. One can conclude that for a less massive planetary compan- 
ion the line profile distortions are weaker. The same conclusion 
can be obtained by the growing influence of the increasing plan- 
etary mass on the disk eccentricity, see Fig. HI 

A natural way to study the effect of the mass of the cen- 
tral star on the line profile distortions would be that the plan- 
etary mass is kept fixed, and the stellar mass is increased. On 
the other hand, we recall that dimensionless units were used in 
hydrodynamical simulations thus the mass of the planet is ex- 
pressed in stellar mass units. Thus changing either the planet 
mass or the stellar mass is equivalent to changing the ratio be- 
tween the planetary and the stellar mass. Therefore, one could 
conclude from the previous results that the larger the stellar 
mass, the stronger the line profile distortions. The situation is 
however a bit more complicated; if we change, for instance, the 
stellar mass, the luminosity of the star changes, which has an 
influence on the temperature distribution in the disk atmosphere 
and interior as well. The effect of the stellar mass can be there- 
fore studied if the planetary-to-stellar mass ratio is kept fixed. 



while the mass along with the surface temperature and the stel- 
lar radius is changed in the synthetic spectral model. Due to 
the effect of the increased/decreased flux of stellar irradiation 
in case of the larger/smaller stellar mass, the outer boundary of 
the CO emitting region will be moved farther/closer to the star, 
and increased/decreased in size. To investigate the dependence 
of the observable line profile distortions on the stellar mass, we 
have recalculated the above presented line profiles for a 0.5 Mq 
(models #l-#4) and a 1.5 M© (models #9-#12) central stars. Part 
of the results are shown in Fig. ^c) (for model #4) and Fig. 
jSld) (for model #12), respectively. Comparing the line profiles to 
those presented in Fig.lSla) (for model #8), it is evident that the 
profiles are considerable narrower for 0.5 Mq and broader for 
1 .5 Mq models, due to the change in the Keplerian angular ve- 
locity {Q.k(R) - (Gnit/R^y^^). Moreover, the line-to-continuum 
ratio at the maximum is slightly decreased and increased about 
a same amount for disks with 0.5 and 1.5 Mq star The change 
in Hne-to-continuum ratio can be easily explained: for a smaller 
mass star, the disk temperature is also decreased due to the de- 
crease in stellar luminosity in a way that the ratio of the CO line 
to the disk plus star continuum is decreased too. The opposite 
is true for larger stellar masses. Finally, we can conclude that 
for a given planet-to-star mass ratio the larger the central stellar 
mass, the larger the observable CO line profile distortion caused 
by the giant planet. Here we have to note that although the obser- 
vational data on T Tauri stars ( |Akeson et al. 200 5) show that the 
size of the disk inner cavity is increasing with increasing stellar 
luminosity (which can be the consequence of the increased dust 
sublimation radius, see Eq. Q), resulting in smaller continuum 
and stronger relative line strength, in our cases the size of inner 
cavity stayed fixed in order to make models comparable. 

Analyzing Fig. IHlb) we conclude that the line profile dis- 
tortions in model #5 (mpi = 1 Mj, w, - 1 Mq, 1 AU orbital 
distance) are so small that the giant planet signal below 1 Mj is 
strongly suppressed. Note that a disk with a 1 Mj mass planet 
embedded into it will not get into the eccentricity state at all, 
see Fig. HI Planets orbiting the host star closer than 1 AU, how- 
ever, can induce stronger distortions, as shown in the following 
section. 



4.4. Orbital distance of the planet 

In this section we investigate the effect of the orbital dis- 
tance of the planet to the strength of the spectral line distor- 
tion. One would initially expect that the spectral line distor- 
tions are stronger/weaker as the orbital distance of the planet de- 
creases/increases. We recalculated the hydrodynamical simula- 
tions in "tight" and "wide" models, where the planets orbit at 0.5 
and 2 AU distances from the central star, respectively. In good 
agreement with the expectations, we found that the planet sig- 
nature in the spectral line was strongly suppressed in the "wide" 
versions of models #l-#8, where m« = 0.5 Mq and nit - I Mq. 
The reason for this is obvious: in wide systems the planet is 
orbiting in cold regions (the atmosphere temperature is varied 
between 200-350 K at 2AU, depending on the stellar mass), 
where the atmospheric CO emission measured to the continuum 
is weak. As can be seen in Fig.lSle), the "wide" version of model 
#8 (»jpi = 8Mj, m» = 1 Mq, 2 AU orbital distance) shows con- 
siderably weaker planet signatures in the CO line profile than 
the original 1 AU model, see Fig.lHla) for comparison. The dis- 
tortion patterns appear at lower velocities owing to lower orbital 
velocities in 0.5 M^ models. On the other hand, we found about 
the same strength of the planet signal in the "wide" version of 
model #12 {m^\ - 12 Mj, m» = 1.5 Mq, 2 AU orbital distance) 
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mass, and indeed this is confirmed by observations ( [Akeson et al.| 
2005 y As we mentioned in Sect. 3.2, the size of the inner cavity 



^ (AU) 



Fig. 7. Azimuthally averaged eccentricities in the "wide" (dotted 
black) and "tight" (solid black) version of model #8 (wpi = 8 Mj, 
m, - \ Mq) and the atmospheric temperature (red). It is appar- 
ent that even though the eccentricity maximum is smaller in the 
"wide" system, the disk is more eccentric on average than in 
the "tight" model. At the same time the planet signal becomes 
stronger in the "tight" (see Fig.lStf)) than in the "wide" (see Fig. 
[8]^e)) model, because the disk atmosphere is considerably hotter 
at 0.5 AU than 2 AU, where the disk is dynamically perturbed by 
the planet. 



as in model #8. This is expected, since in the latter case the lu- 
minosity of the host star is high enough to produce sufficient CO 
excitation even at 2 AU. 

The resulting line profile calculated in the "tight" version of 
the model #8 (»v = 8 Mj, m» = 1 Mq, 0.5 AU orbital distance) 
is shown in Fig.lSjf). Because the giant planet orbits closer to the 
host star, where the disk atmosphere is heated to higher temper- 
atures (the atmosphere temperature varies between 450-700 K at 
0.5 AU, depending on the stellar mass), the emission emerging 
from the perturbed regions is stronger than in the original 1 AU 
models. It is also evident that the distorted peaks are slightly 
shifted to higher velocities compared to the original 1 AU mod- 
els, because the orbital velocity of the gas parcels in dynamically 
perturbed orbits is higher in "tight" models. The disk is more ec- 
centric on average in the latter case (Fig. ITli, but the strength of 
the line profile distortion is increasing as the planetary orbital 
distance decreases. 

In order to demonstrate that smaller mass planets can also 
be detected in closer orbits we have calculated the line profile 
in a "tight" version of model #2, m^i = 1.5 Mj, m» - 0.5 Mq, 
0.5 AU orbital distance. We find that one could observe roughly 
the same strength line profile distortions as in model #8 (wpi - 
8 Mj, m^-l Mq, 1 AU orbital distance). 



has a significant influence on the overall line-to-continuum ratio. 
In disks dynamically perturbed by planets orbiting at 1 AU and 
with a significant amount of CO lying inside 0.2 AU, the ma- 
jority of CO flux is emerging from the unperturbed innermost 
regions, which eventually could smear out the planet signatures. 
On the contrary, if the size of the inner cavity is larger, the con- 
tribution of the perturbed regions to the total CO flux becomes 
stronger In Fig. [8 g) we present the line profiles obtained in 
model #8, in which the inner cavity is increased to 0.4 AU in 
size, while the 8 Mj mass planet was still orbiting at 1 AU. As 
can be seen, the overall line-to-continuum ratio is decreased due 
to the absence of hot gas, but the level of the asymmetric pattern 
profile is more apparent than in models with an inner cavity ex- 
tending only to 0.2 AU, see Fig.[8|a) for comparison. For larger 
inner cavities the continuum level of the disk interior and inner 
rim are decreased too, resulting in a moderate weakening of the 
line-to-continuum ratio and slight strengthening of planet signal. 
In this case the detection possibility of a planet orbiting at larger 
distances (> 2 AU) is growing. 



4.6. Disk geometry 

It is well known that many TTauri sources present flatter than 
AF,x ~ A^'^'^ SEDs. One attempt to expl ain fliis SEP flatten- 
ing was to assume that the disk is flaring (Kenyon & Hartmann| 
19871, resulting in a geometrically thicker disk atmosphere 
H(R) ~ R'>''^\ where y is the flaring index. In the approach of 
Chiang & Goldreich ( 1997 1 the flaring index is y = 2/7. In order 
to investigate the effect of disk geometry on the strength of the 
line profile distortions we recalculated several models using flat 
thermal disk approximation (the hydrodynamical inputs were the 
same). An example of the flat version of model #8 (mpj - 8 Mj, 
m» - 1 Mq, 1 AU orbital distance) is shown in Fig. Imh). Note 
that in order to have the same peak line-to-continuum ratio in 
flat as in flared models, the amount of emitting CO is arbitrary 
increased by decreasing the dust-to-gas ratio to 1.2 X lO""* and 
keeping the disk mass unchanged. The reason for this is that al- 
though the atmospheric temperature does not differ in flared and 
flat disk models (see the atmospheric temperature given by Eq. 
( |A.5| l), the atmospheric surface density of CO is substantially 
smaller in flat than flared disk models (see the dependence of 
surface density given by Eq. (D.2i on the grazing angle, given 
by the first term of Eq. 



( |A.2[ ) in flat disk models). It is evident 
that the total line-to-continuum ratio is higher in the flat than in 
the flared disk model, but the line width at the foot of the line 
profile does not change. According to our calculations the line 
profile distortions in flat models compared to the flared one are 
less significant in "wide" systems where planets orbit at larger 
distances. This can be explained by the effect of flaring getting 
stronger with increasing distance to the host star. 



4.5. Inner cavity size 

Now we turn our attention to the question under which circum- 
stances planets can be detected orbiting farther away from the 
star (farther than 2AU for instance). In Sect. 3.2 we argued 
that there is a region around the star called inner cavity from 
which the CO is already depleted. In our simulations presented 
so far, we assumed that the radius of the inner cavity is fixed at 
0.2 AU irrespective of the stellar luminosity. On the other hand, 
the size of the inner cavity is presumably dependent on the stellar 



4.7. Observational considerations 

To observe the planet-induced distortions of the CO line pro- 
files that we simulated, a spectroscopic facility is required that 
meets two basic demands: it needs sufficient spectral resolution 
to properly resolve the sub-structure in the line profile and suf- 
ficient sensitivity to bring out the low-contrast line profile dis- 
tortions. In this section, we investigate the observability of the 
modeled effect in the context of contemporary facilities, and give 
a brief outlook into the ELT era. We present a quantitative ex- 
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Fig. 8. Distorted CO line profiles emerging in different models of planet-bearing disks. The disk inclinations were 20°, 40° and 60° 
shown in blue, green and red colors. For comparison, we display the symmetric double-peaked line profiles (dashed lines) emerging 
from a planet-free disk as well. Line profiles shown in panel (a) and (b) show the dependence of the distortion strength on the 
planetary mass. In panels (c) and (d) the effect of stellar mass, in panels (e) and (f) the effect of the orbital distance is presented. The 
effect of the inner cavity size on the strength of the line profile distortions is presented in panel (g). The impact of the disk geometry 
is shown in panel (h), where the line profiles emerging in model #8 were calculated in flat-disk geometry. Note that the line profiles 
emerging in distinct models is calculated at orbital phases where the line profiles asymmetry is the most prominent, i.e., the orbiting 
planet is not always at the same position angle in each figure. 
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Fig. 9. Permanent asymmetry and variations in the distorted CO ro-vibrational line profiles in a disk with an 8 Mj planet orbiting 
an 1 Mq star at 1 AU (model #8). The subsequent epochs go across and down. The position angle (PA) of the giant planet is shown 
in each frame. The arrows point to the variable component moving in between ^ ±25km/s. The normalized line flux is shown 
with the same scale as in Figj8|a). The snapshots were calculated during one planetary orbit, thus about 2.5 weeks elapsed between 
the frames. The disk inclinations were 20°, 40° and 60°, shown in blue, green and red colors. For comparison, we displayed the 
symmetric double-peaked line profiles (dashed lines) emerging from a planet-free disk as well. As can be seen, the line profiles have 
permanent asymmetry and show a significant change in their shapes within weeks. 



Zs. Regaly et al.: Detectability of giant planets in disks 



13 



m^-lMg, mp,-8Afj, <3p,=lAU 




Fig. 10. Effect of instrumental noise and atmospheric absorption 
on the observability of line profile distortions (top) and the Earth 
atmospheric transmission at 4.7545yum (bottom). The line pro- 
files were calculated in model #8, shown already in Fig.|8ja), but 
superimposed with an artificial noise expected for one hour of 
integration on M^ag = 6 brightness source with VLT/CRIRES. 



ample for the Cryogenic High Resolution Echelle Spectrograph 
(CRIRES) at the VLT (KaeufI et al. 2004), which is currently 
the most powerful high-resolution spectrograph covering the M- 
band around 4.8 yum, capable of observing CO ro-vibrational 
spectra. 

The main limiting factor for ground-based thermal infrared 
observations is the Earth atmosphere, which strongly absorbs the 
radiation from astronomical sources and causes a high thermal 
background. Particularly strong telluric absorption and emission 
is present at the wavelengths of low-excitation CO transitions. 
Within ~ lOkm/s from the line center the telluric transmission 
is so low and the background emission so high that these spec- 
tral regions are not accessible. However, depending on the loca- 
tion of the source with respect to the ecliptic, the CO features 
are Doppler shifted by up to ±30km/s due to the Earth's orbital 
motion, making the whole CO line accessible over the course 
of several months, in principle. Note though that we can never 
measure the entire CO line profile simultaneously, there will al- 
ways be an approximately 20km/s wide "gap" in our spectral 
coverage at any observing epoch. 

In order to test whether we would be able to detect the planet- 
induced permanent profile asymmetries with currently available 
instrumentation, we have simulated a VLT/CRIRES observation 




Fig. 11. Goodness-of-fits of distorted line profiles contaminated 
with artificial noise obtained in planet-free (dashed lines) and 
giant planet-bearing disk (solid lines) models versus the integra- 
tion time for 20°, 40° and 60° inclinations, represented by blue, 
green and red colors, respectively. It is visible that for a reason- 
able signal-to-noise ratio, the planet-free models provide always 
fits with less confidence than the giant planet models, i.e the 
X^ is always larger Moreover, the smaller the inclination angle, 
the larger the required minimal integration time (fint > 3187 sec, 
1909 sec and 2271 sec for the 20°, 40° and 60° incHnations, re- 
spectively) to distinguish the planet-free and planet-bearing disk 
models. 

of model #8. We have assumed the source to have a brightness 
of M = 6 mag, and calculated the achieved SNR as a function 
of wavelength using the VLT/CRIRES exposure-time calcula- 
tor] This calculation takes into account the telluric absorption 
and emission, the quality of the adaptive optics correction, the 
system throughput and detector characteristics. We note that our 
simulation only considers the signal-to-noise ratio, and assumes 
that systematic effects due to, e.g., time-variable telluric absorp- 
tion, can be well calibrated. The achieved SNR can be scaled 
according to 



SNR cc 10' 



,0.4M 



(8) 



where M denotes the M-band magnitude and Tint the integration 
time. In Fig.[TO]we show the emerging V=1-0P(10) line profiles 
for model #8 (upper figure) and the atmospheric transmission 
(lower figure) for an assumed integration time of 1 hour on a 
source of M = 6 mag. The model spectrum has been convolved 
with the 3 km/s CRIRES instrumental resolution. We have as- 
sumed zero relative velocity between the source and the Earth at 
the time of the observation. The effects of the telluric CO line are 
clearly seen around the line center, where no meaningful mea- 
surement can be obtained. The line wings, where the effects of 
the planet are most prominent, can be well measured. 

In order to quantify the significance with which a plane- 
tary signal (the most prominent permanent asymmetric pattern) 
may be detected, we have tried to fit the modeled observation 
with both planet-free and planet-bearing models. In Fig. 11 we 
show the goodness-of-fits (reduced ;t'^) as a function of integra- 
tion time for both cases (with dashed line for planet-free and 



"* We used version 3.2.8 of the ETC, assumed a water vapor column 
of 2.3 mm, an airmass of 1.4, a seeing of I'.'O, a slit width of 0'.'2, and 



adaptive optics correction using a guide star of R ■■ 
tral type of MOV. 



10.0 mag and spec- 
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solid line for planet-bearing disk models) for 20°, 40°, and 60° 
inclination angles. Above a certain integration time (presented 
by vertical lines) the planet-bearing models yield a more than 
50% better reduced x^ fit to the simulated data than those of the 
planet- free models. The asymmetric pattern induced by an em- 
bedded planet increase and thus become easier to detect with 
disk inclination. 

Because the fundamental band ro-vibration R(0), R(l), 
R(12), R(13), R(19)-(21), R(25), R(26) and R(30)-(33), P(6-12) 
and P(22-32) CO lines are not blended by the higher vibrational 
linea^ the line profiles can be averaged, using appropriate scal- 
ing, resulting in an increased signal-to-noise ratio. Note that the 
unblended lines with higher 7 > 25 or lower J < 2 rotational 
quantum numbers are too weak for averaging: using them would 
introduce only additional noise to the averaged profiles. Also 
note that, since the higher vibrational levels are not excited at 
all in disk atmospheres similar to that used in our models, the 
blended lines can also be used for averaging, resulting in an even 
more increased signal-to-noise ratio. 

We conclude that in the confines of our model, the perma- 
nent asymmetric signals produced by an 8 Mj giant planet orbit- 
ing at ~ 1 AU embedded in the disk of a T Tauri star of bright- 
ness M - 6 mag would be detectable with VLT/CRIRES with an 
integration time of 1 hour Moreover, because the variable com- 
ponent at the top of the permanent asymmetry has a width of 
~ 5 - lOkm/s, exceeding the CRIRES ~ 3 km/s resolution, and 
~ 10% magnitude, there is a certain possibility to strengthen its 
planetary origin. 

With the next generation of giant, ground-base telescopes 
our sensitivity to low-contrast features in infrared spectra will 
increase greatly. We have used the ESO exposure-time calcula- 
tor for the ELT (version 2.14) to simulate an observation with 
the proposed first generation infrared instrument Mid-infrared 
E-ELT Imager and Spectrograph (METIS) (B randl et al.|2008] l. 
The exposure-time calculator necessarily incorporates planned 
telescope and instrument specifications, because this facility is 
yet to be built. We find that ELT/METIS yields an increase in 
sensitivity of a factor of ~30 compared to VLT/CRIRES. Thus 
also sources that are substantially fainter or have weaker spectral 
features than those modeled here in the context of VLT/CRIRES 
are within reach of the observational facilities of the next decade. 



5. Conclusions 

We have shown that our semi-analytical 2D double-layer thermal 
disk model (in which the disk layers are heated mainly by stellar 
irradiation, and the emission of the disk inner rim is accounted 
for) is capable to reproduce the double-peaked Keplerian line 
profiles in CO fundamental ro-vibrational band in young (about 
2.5 Myr old) T Tauri-type protoplanetary disks assuming canon- 
ical dust and gas properties. However, note that several T Tauri 
stars show centrally peaked symmetric profiles instead, which 
can be attributed to an enhanced CO flux formed in distant re- 
gions (up to 10 AU), presumably above the disk atmosphere. 
This can be explained by non-LTE heating, such as UV fluo- 
rescence (jKrotkov et al.||1980), or dust gas temperature decou- 
pUng (j^assgold et al. 2004; Kamp & Duflemond 2004), result- 
ing in highly excited CO emissions. These effects were not taken 
into account in our calculations. Thus the gas temperature at 
disk atmosphere provided by our thermal model is presumably 



^ Assuming that the abundance of '-C"'0 isotopologues are low, thus 
their fundamental V = 1 — > contribution to lines is negligible. 



somewhat lower than the realistic one. Our predictions, how- 
ever, demonstrate a "the worst case" because the planet signal 
is getting stronger with increasing atmospheric temperature. It 
is revealed that significant line profile distortions appear in the 
CO fundamental ro-vibrational band caused by a giant planet. 
As we have demonstrated, a giant planet with a mass of at least 
1 Mj can perturb the disk sufficiently to be detected by perma- 
nent asymmetry in line profiles and its short timescale variations. 
The permanent asymmetry in line profiles can be interpreted by 
the disk being in an eccentric state only in the gap, but the short 
timescale variability is certainly connected to the local dynami- 
cal perturbations of the orbiting planet. Depending on the signal- 
to-noise ratio and the resolution of the spectra, even lower-mass 
giant planets can be detected revolving on close orbits in disks 
observed with high-inclination angles and with large size inner 
cavities. The summary of our findings is: 

1 . The dynamical perturbation induced by a giant planet signif- 
icantly distorts the CO ro-vibrational fundamental emission 
line profiles. 

2. The line profile distortions are more apparent for larger in- 
clination angles. 

3. The main component of line profile distortions is a perma- 
nent asymmetry. Its position in the line profile depends on 
the stellar mass and orbital distance of the planet. The per- 
manent asymmetry shape depends on the orientation of the 
gas elliptic orbit with respect to the line of sight. 

4. The line profiles are changing with time, correlated to the 
orbital phase of the giant planet. The timescale of the vari- 
ation is on the order of weeks, depending on the orbital pe- 
riod of the planet, i.e, the orbital distance and the mass of the 
host star. The magnitude of variable component is ~ 10%, 
its width is ~ 10 km/s, depending on the inclination angle. 

5. The planet signal becomes stronger with increasing plane- 
tary mass. 

6. The planet signal strengthens/weakens with increas- 
ing/decreasing mass of the host star, neglecting the effect of 
the stellar luminosity on the size of the inner cavity 

7. The planet signal strengthens/weakens with decreas- 
ing/increasing orbital distance of the planet. 

8. The size of the inner cavity has a strong influence on the gi- 
ant planet observability. If the size of the inner cavity is sub- 
stantially smaller than 0.2 AU, we can detect a giant planet 
only at close (about 0.5 AU) orbit. On the contrary, if the 
cavity is larger (0.4 AU), a distant planet (in 2 AU orbital dis- 
tance) can also be detected. 

9. The influence of disk geometry on the planet signal is mod- 
est. With a decreasing flaring index the strength of planet sig- 
nal does not change although the line profiles are strength- 
ened, except in "wide" flat models where the planet signal 
substantially suppressed. 

10. The lowest mass of giant planet orbiting at 1 AU that still can 
be detected is 0.5 Mj for a 0.5 Mq star in the confines of our 
models. 

11. By high-resolution near-IR spectroscopic monitoring with 
VLT/CRIRES, giant planets with > 1 Mj mass orbiting 
within 0.2-3 AU in a young disks may be detectable, if other 
phenomena do not confuse, mimic, or obscure the planet sig- 
nature. The level of confidence is growing with increasing 
inclination angle. 

In the light of our findings, we propose to observe T Tauri 
stars that shows an appreciable amount of CO in emission 
with high-resolution near-IR spectrograph to search for giant 
planet companions. The asymmetric and time-varying double- 
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peaked line profiles can be explained by strongly non-circularly 
Keplerian gas flow in the disk caused by the giant planets, as 
was presented. Although higher V > 2 vibrational states of CO 
are not excited in our models, overlaps can occur, resulting in 
asymmetric line profile. To avoid this, one should consider only 
non-blended lines for profile averaging such as P(6)-P(12), and 
R(0), R(l), R(12), R(13), R(19)-(21), R(25), R(26) and R(30)- 
(33). Another possible source of permanent line profile distor- 
tions could be a large companion deforming the circumprimary 
disk to fully eccentric in mid-separation binaries (e.g. Kley et al. 
( 2008| l; Regaly et al. in prep.). Last but not least, the distorted 
PSF, caused by tracking errors of the telescope or unstable ac- 
tive optics during an exposure, can induce artificial signals. A 
real signal from a bipolar structure will change sign when ob- 
served using antiparallel slit position angles whereas the artifi- 
cial signal will not. Thus artificial signatures can be successfully 
identified with the comparison of the spectra of disks obtained at 
two antiparallel slit position angles, e.g. 0° and 180° (Brannigan 
2006| l. By detecting short time-scale (weeks to months) 



et al. 



variations of the asymmetric profiles with similar patterns as pre- 
sented in Fig. l9j we could infer the presence of a giant planet in 
formation. 



The proposed instrument METIS for the ELT (Brandl et al. 
2008 I, providing high-resolution (AjhA - IQr') infrared spec- 
troscopy with an increase in sensitivity of a factor of ~ 30 to 
that of VLT/CRIRES will give the opportunity to even decrease 
more the mass and brightness limits of possible giant planets. If 
we revealed for instance the distance of the observed planets to 
the star as a function of the age of the star-disk system, we would 
be able to constrain planetary migration models. Moreover, de- 
termining the youngest star host planet-bearing disks could also 
provide constraints on planet-formation models. 

Although our assumptions are reasonable from the point of 
view of finding a simple model of the signal of dynamical per- 
turbations, it can be developed further of course. For example, 
if the gap is cleared by the planet to a fraction of the surround- 
ing density, the disk is not directly irradiated, but is shadowed by 
the inner wall of the gap ( Varniere et al. 2006 ). This mild heating 
causes a temperature drop at the inner gap wall, which results in 
a depressed CO flux. On the other hand, the outer edge of the 
gap is directly illuminated by the star, resulting in a slight tem- 
perature increase, which produces stronger CO emission from 
the regions where the eccentricity reaches its maximum (see 
e.g. FigHl. According to our results the gap is strongly non- 
axisymmetric, thus 3D radiative transfer using 3D density struc- 
ture needed to calculate the real temperature profile in perturbed 
disks. Moreover, as the gap is depleted in large grains because 
of the trapping of dust grains larger than ~ 10 /^m ( |Rice et aT] 



2006 1, the CO flux normalized to the continuum is strengthened. 



The gas and dust can be thermally uncoupled (see e.g., Glassgold 
[eTai. (2004); Kamp & DuUemon d (2004); Woit ke et al.| ( |2009| l) 
in the tenuous gap, which may result in CO being overheated 
to the dust. These effects could cause additional permanent line 
profile distortions due to the elliptic shape of the gap seen in den- 
sity distribution (Fig.[3]l. If we consider an excess emission orig- 
inating from the planetary accretion flow suggested by Clarke 
|& Armita ge ( 2003 ), or the density waves in disk surface caused 
by the orbiting planet, it would cause a varying planet signal. To 
these caveats, one might also consider the disk turbulence (e.g., 
due to MRI effects), which may add random asymmetries and 
"blobbiness" to the line profiles. We will investigate these ef- 
fects in a forthcoming paper using a more sophisticated thermal 
disk model. 



CO exists inside the dust sublimation radius in some T Tauri 
stars ( |Carr||2007| l or Herbig Ae/Be stars dBrittain et"aL]|2009| l, 
where the disk is optically thin in the absence of dust. One pos- 
sibility to excite CO in the optically thin inner disk close to the 
stellar surface is the UV fluorescence, which is far from LTE 
(Krotkov et al. 1980[ l applied in our calculations. Depending 
on the efficiency of X-ray heating the gas temperature can be 
as high as ~ 4000 - 5000 K. Somewhat below this tempera- 
ture, where CO still exists, the high-excitation vibrational lev- 
els of CO such as y > 2 are populated, resulting in significant 
y = 2^ l,y = 3^2, etc. ro- vibrational lines. According to 



Skrutskie et al. ( 1990), the inner cavity in transitional disk may 



be formed by a close giant planet. Indeed, Lubow & D'Angelo 



(2006 1 found that an 1 - 5Mj mass planet significantly lowers 
the accretion of dust and gas to 10% -25% of the accretion rate 
outside the gap. We expect that the planetary signal forming in 
the optically thin inner disk is also detectable, which is strongly 
supported by the revealed line profile distortion strengthening 
with decreasing orbital distance of the planet. The UV pumping 
requires considerable UV flux from the central star, which is a 
characteristic of Herbig Ae/Be stars rather than smaller T Tauri 
stars. Considering that the CO may be excited up to to lOAU 
by UV fluorescence ( |Brittain et al.|2007| [2009, ) in Hebig Ae/Be 
stars, planet orbiting at larger distances than 2 AU may be also 
detectable by CO ro-vibrational line profile distortions. 

Spectroscopy of T Tauri stars shows emission of molecules 
such as H2O, OH, HCN, C2H2, and CO2, as well. Nevertheless 
CO is more abundant than these molecules by ~ 10 times, only 
H2O could reach the abundance of CO, predicted by recent mod- 
els that calculate the vertical chemical structure of the gas in disk 
atmosphere (e.g. [Glassgold et"ar ( |2004| l, |Kamp & Dullemond] 
( |2004| i, and Woitke et al. (|2009f ). Indeed in some cases, like 
AATau (|Carr & Najita|2008 |), and AS 205 A and DR Tau ( (Sal>¥ 



|et al.|2008| l, the rotational transitions of H2O dominate the mid 
infrared (10 - 20;um) spectra, suggesting that H2O is abundant 
in disk atmospheres. The strong water emission could be the 
consequence of turbulent mixing that carries molecules from the 
disk midplane, where they are abundant, to the disk atmosphere 
(Carr & Najita 2008), or the effects of an enhanced mechani- 
cal heating of the atmosphere ( [Glassgold et al. 2009). While the 
H2O ro-vibrational lines probe the inner regions of disk out to 
radii > 2 AU, the rotational lines are produced between radii 10- 
100 AU (Meijerink et al. 2008) 1. The rotational and ro-vibrational 
line profiles of water are also subject to distortions caused by 
the dynamical perturbations of a giant planet. Thus it is reason- 
able to search for signature of giant planets in the high-resolution 
spectra of H2O in the mid-infrared band (for planets orbiting int 
the inner disk), and rotational spectra of H2O in the far-infrared 
band (for planets orbiting in the outer disk), as well. Because 
H2O is heated by stellar X-rays and sub-thermally populated 
beyond 0.3 AU, X-ray-heating and non-LTE level population 



treatment is needed to calculate water lines (e.g., Meijerink et al. 
([2008'); Kam p et al.| ( |20T0l l) 



Finally, we point out some noticeable resemblance to our 
findings revealed in observed CO line profiles. The V836Tau 
transitional disk shows strongly distorted 4.7 jum CO ro- 
vibrational line profile presented by Najita et al.[ ([2008 ). The 
averaged CO line profiles shows very similar features to our 
results. According to jNajita et al. (2008i the possibility of a 



massive planet is also restricted by current limits on the stel- 
lar radial velocity of V836 Tau, which constrain the mass of a 
companion within 0.4 - 1 AU to 5 - 10 Mj. Because the disk 
extents to a limited range of radii (~ 0.05 - 0.4 AU), the dynam- 
ical perturbation of a close orbiting planet with a mass smaller 
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than 5 Mj presumably could cause the observed line profile dis- 
tortions. Note that radial velocity measurements indicate that 
no planet is present larger than 1 - 2 Mj closer than ~ 0.5 AU 
( jPrato et al.||2008| . The transitional disks that encompass the 
young sources SR21 and HD 135344 B have been observed by 
VLT/CRIRES UPontoppidan et al.|2008| and they already show 
clear line profile asymmetries in CO fundamental ro-vibrational 
band. Our scenario of a planet is a reasonable explanation of 
the asymmetry. The presence of a planet was also proposed by 
Grady et al.| ( |2009) for HD 135344 B and' Eisner et aL] ( [2009l ) for 
SR21, based on SED and visibility data, respectively. EXLupi 
(prototype of EXor type young variable stars) showed remark- 
able line profile variations during its 2008 outburst (Goto et al. 
2010, in prep.). Goto et al. revealed that the CO emission has a 
spatially multiple origin. The quiescent component forms in the 
outer optically thick disk, while the outburst component in the 
inner optically thin gas disk, where the higher vibrational lev- 
els of CO are presumably excited by UV pumping. These high- 
excitation y = 2 — > 1 and V = 3 — > 2 lines are double peaked, 
strongly asymmetric and variable on a short (weekly) time scale. 
Thus for example, by measuring short timescale (on the order of 
weeks or less) periodic variations in the CO ro-vibrational spec- 
tra of these sources, by mid-lR monitoring observations would 
eventually indirectly detect an embedded Jupiter-like planet in 
birth. 
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Note that for a non-flared, i.e flat disk the grazing angle is de- 
fined by only the first term of Eq. ( A.2 1. 

By definition the optical depth of the disk atmosphere along 
the incident stellar irradiation is ry = 1 in the optical wave- 
lengths. The optical depth of the atmosphere at near-IR heated 
to ratm,in(^) perpendicular to the disk plane is Ty6{R)£atm(R), 
where eatm(^) - iT,itm,mi^)/T*y^ is the dust emissivity, and /3 
is the power law index of the absorption coefficient of the dust, 
see Appendix [C] for details. The total flux emitted inward and 
upward by the optical thin atmosphere can be given by 



F,UR) = 26{R)£,,^(R)crT,,^,„(R)\ 



(A.4) 



In LTE the absorbed flux of atmosphere equals to the emitted 
one (Ft{R) - Fatm,irr(^)), thus the disk atmosphere is heated to 
the temperature 



^ atm,in-(") — 



j,l/(4+/3) ^ ^2/(4+/J) 



m 



T, 



(A.5) 



by stellar irradiation. 

Assuming that the disk interior with a temperature Tinii„(R) 
is optically thick to its radiation everywhere in our computa- 
tional domainr] the flux F[„t,[„(R) emitted by the disk interior 

is 



■Fint,in-(^) - 0'7'int,in-(^) 



(A.6) 



Considering that the absorbed stellar flux will be re-emitted by 
the disk atmosphere upward and downward only half of this 
flux heats the interior, i.e 0.5Fatm(^) = O.SF^R) = Fi„t,i,T(^). 
Accordingly, the temperature of the disk interior set by the stel- 
lar irradiation can be given by 



Appendix A: Temperature distribution in the disic 

For the sake of completeness we present here the derivation of 
the dust temperature distributions in the disk interior and atmo- 
sphere according to 'Chiang & Goldreich (1997). Below we use 
the local thermodynamical equilibrium (LTE) assumption every- 
where. The stellar flux in the disk atmosphere at a distance R 
from the stellar surface is 



«.).^(^f-. 



(A.l) 



where 6{R) is the grazing angle of the incident irradiation, Rt and 
Tt are the stellar radius and surface temperature, respectively, 
and cr is the Stefan-Boltzmann constant. Here we assumed that 
only half of the stellar surface is visible from a given point at 
a distance R from the star in the disk atmosphere. According to 
|Chiang & Goldreich ( 1997 1 the grazing angle can be given by 



--^li^yMii) 



-2/7 



(A.2) 



in a flared disk assumed to be in hydrostatic equilibrium. In Eq. 
( A.2 1 Tg is the gravitational temperature at which the thermal en- 



ergy of gas parcel balances the gravitational energy at the stellar 
surface, which can be given by 



G Merrill 
kR. ' 



(A.3) 



where ;«h is the mean molecular weight of the gas, G and k is 
the gravitational constant and Boltzmann constant, respectively. 



r....^{^m'r.. 



(A.7) 



Owing to the accretion process a significant amount of grav- 
itational potential energy has to be dissipated by a viscous pro- 
cesses. According to |Lynden-Bell & P ringle ( 1974), in a steady 
state disl|jwith a constant accretion rate M, the flux F^^cdR) re- 
leased by the disk mid-plane due to the change of potential en- 
ergy can be given by 



;.„(« = 5£^„ 



-(I)'" 



R 



(A.8) 



In this way the disk mid-plane is heated to the temperature 



T,cc(R) 



3GM,M 



(I)" 



im-|i/4 



R-''\ 



(A.9) 



Here we have to note that only half of the total accretion power 
is involved in Eq. ( |A.8[ ), the remaining is stored in the kinetic en- 
ergy of orbiting gas. This energy should be radiated away by the 
disk boundary layer (Popham et al. 1993 , 1995) or by the accret- 
ing material in the funnel flow formed along the magnetic field 
lines in magnetospheric accretion model, see Hartmann et al. 
( 1994 1; Bouvier et al. (20071. Because the disk boundary layer 
and the funnel flows are confined into such small volumes that 



* As mentioned before, in the computational domain, the disk interior 
remains optically thick to its own radiation. 

' The disk being in steady state means that the density distribution 
corresponding to the surface density distribution does not considerably 
change in time. 
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the gas temperature reaches about lOOOOK emitting in the UV 
band, its radiation is not taken into account. 

Now let us take into account the heating due to viscous dis- 
sipation in disk interior using the superposition principle, i.e the 
radiation at a given location is the sum of the radiations cor- 
responding to different heating sources. To determine the disk 
interior temperature Ti^j^cdR) caused by the viscous dissipation 
first assume that the flux emitted by the optical thick disk interior 
is 



■f'int,acc(^) - 0'7'int,acc(^) • 



(A. 10) 



Because the disk mid-plane radiates the accretion flux into two 
directions 0.5Facc(^) = ^int,acc(^), the disk interior is heated to 



Tnn,..c(R) = (l/2)"^r,ecW 



(A. 11) 



by accretion. Taking into account the heating of the disk interior 
by irradiation of atmosphere and viscous dissipation together 
(simply summing the radiation fluxes), the resulting temperature 
of disk interior is 



1/4 



TiatiR) - [Tiai^iniR) + Tint.accC^) ) 



(A.12) 



Appendix B: Emission from tlie disl< inner rim 

To determine the temperature profile of the disk interior close to 
the disk inner edge, where the simple double layer assumption 
cannot be applied, we first assume that the rim is a perfect verti- 
cal wall. Any given optically thick vertical slab with a thickness 
of dR at R + dR distance from the central star is irradiated by 
the neighboring hotter slab located at R. Assuming that half of 
its radiation is received by the one ai R + dR, the emitted and 
irradiated fluxes are in equilibrium in LTE, i.e. 



1 



(TT,URf = cr€iR)T,UR + dR)\ 



(B.l) 



where €(R) - (Tnm{R + dR)/Tii„^{R)'f is the emissivity of the 
slab at R + dR, see in Appendix |C] Note that here we neglect 
the irradiation of the neighboring slab atR + 2dR with a lower 
temperature than the slab at R + dR. To calculate T^m{R), we first 
approximate T^^^iR + dR) with Tiitn{R) + T'. {R)dR. This results 
in an ODE for T{R), which has the solution 



rrim(^0)exp 



l/(4+«\ 



(R - Rq) 



(B.2) 



where Rq is the radius of the disk inner edge. To take into con- 
sideration the additional irradiation of the rim at the opposite 
side, we assume that the temperature at the innermost slab is 
Tdm(Ro) - qTi„ii„{Ro), where ^ ^ 1.2 according to results of 
our simulations done by 2D RADMC (Dullemond & Dominik 
|2004) . Thus incorporating the additional heating by the disk rim, 
the disk interior temperature, previously given by Eq. (A.12i, 
can be given by 

TnniR) = {TintMRt + rint,accW' + T'rimW')'^^ - (B.3) 

where the superposition rule is applied. 



Appendix C: Dust emissivity 

In this section the derivation of dust emissivity at a specific tem- 
perature is given, assuming optically thick environment. Let us 
define the dust emissivity at a temperature Tdust heated by a 
blackbody of the temperature Tin- as the ratio of the Rosseland 
mean opacities 

<'f(7'dust))R 



The Rosseland mean opacity is defined by 



{<T)}^ = 



J^ dB(v, T)/dTdv 
J^ dB(v, T)ldTK{v, T)-^dv 



(C.l) 



(C.2) 



where B{v, T) is the Planck function. The dust used in our model 
generates the k{v) - Koiv/v^'f opacity law, where /3 > 0. Note 
that wee use P - I throughout all our models (Rodmann et al. 
2006). Following the calculations of Stabl er & PaUa| ( |2005| l in 
Appendix G, we assume that the dust absorption coefficient can 
be given by 

K{v,T)^Ka(^\ /, (C.3) 

where k and h are Boltzmann and Planck constants, respectively, 
and X - hv/kT. Supposing that 



dB(v, T)/dTdv = Af{x)dx, 



(C.4) 



where A is an appropriate dimensional constant, and substituting 
this into the definition of Rosseland mean opacity, Eq. (C.2i, 
leads to 

{4T)h^Kj — \ / . 

yhvQJ j x-f^f(x)dx 

Because the quotient of integrals is a pure numb er, th e R ossela nd 
mean opacity is proportional to T^. Applying 
find that the dust emissivity can be given by 



(C.5) 



C.l 



and 



C.5 



we 



fdust — 



T'dii 



(C.6) 



Appendix D: Optical depth of the disk atmosphere 

Assuming that the disk is observed with an inclination angle /, 
the monochromatic optical depth of the disk atmosphere is the 
sum of optical depths of dust and gas along the line of sight, i.e.: 

T(v, R, cp, i) = ^— (^d(v)£dW + %(v, R, 4>, i)^g(R)) , (D. 1) 
cos(0 ^ "J 

where Ki(v) and k„(v,R, (p, i) are the dust and gas opacity at fre- 
quency V, while t.d{R) and Eg(/?) are the dust and gas suiface 
densities in the disk atmosphere. Note that because the dust 
opacity (k^) is taken to be uniform throughout the disk, the op- 
tical depth of the disk atmosphere depends on R and (p via the 
gas opacity characterized by the temperature distribution and the 
line-of-sight component of the orbital velocity of gas parcels. 
By definition the optical depth at the bottom of the disk atmo- 
sphere is unity at optical wavelengths along the stellar irradi- 
ation. Assuming that the mere opacity source at visual wave- 
lengths is the dust, the dust surface-density of disk atmosphere 
is 



^d(R) 



m 

Ky 



(D.2) 
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where /cy is the overall opacity of dust at visual wavelengths. 
Assuming that the gas- and dust-mass ratio to the total mass (Xg 
and Xd, respectively) is constant throughout the disk (i.e. there 
are no vertical or radial variations in the mass ratios), the surface 
density of gas in the disk atmosphere is 



Ah 



(D.3) 



using the dust surface-density Eq. (D.2 1. Thus, the optical depth 
in the disk atmosphere along the Une of sight can be given by 



T(v,R,^,i) 



1 / ..6(R) X,6(R) , „ . .. 
7^ kd(v) + Kg(v,R,(f>,i) 

COS(l) \ A-y Ad A-y 



(D.4) 



where we used Eqs. (D.l|D.3|l 



Appendix E: Monochromatic opacities 

In LTE the monochromatic opacity of the emitting gas at fre- 
quency V with a molecular mass of mco due to transitions be- 
tween states u— >1 can be given by 



K„{v,R,(f>,i) 



1 



1 






K\gu 



X exp 

X <b{v,R,4>,i) 



- exp 



kT,,^(R) 



(E.l) 
(E.2) 



where vo is the fundamental frequency of the transition, 
QiT,ni-n{R)) is the partition sum at the gas temperature ratm(^), 
Aui and ^u ^e the probability of transition (i.e the Einstein A 
coefficient of the given transition) and the statistical weight of 
the upper state, respectively, while c is the the light speed. The 
partition sum can be given by 



Q{T,UR))^Yj^i^^'^ 



Ei 



kT,UR) 



(E.3) 



where ^i and Ei are the statistical weight and energy level of 
the ith excitation state. In Eq. (E.l i the <i>{v, R, (p, i) is the local 



intrinsic line profile originating by the natural thermal and the 
local turbulent broadening acting together If the pressure of gas 
is negligible, the intrinsic line profile can be represented by a 
normalized Gauss function 



^(v,R,<p,i) 



1 



r(R)^/^ 



exp 



Vo + Av(R, (j), i) 



(E.4) 



where criR) is the line width, Av{R, (p, i) is the line center shift 
due to Doppler shift caused by the apparent motion of the gas 
parcels along the line of sight. The line width is determined by 
the natural thermal broadening 



CTthei 



.(R) = - 

c 



2kT,UR) 
mco 



and the local turbulent broadening 



C V '«H 



(E.5) 



(E.6) 



assuming that the speed of turbulent motions is;^' times the local 
sound speed. In Eqs. (|E.5[) and (E.6 1 mco and me are the molec- 



ular masses of H and CO and y is the adiabatic index of the main 



constituent of the gas, i.e that of hydrogen. Considering the ther- 
mal and local turbulent broadening, the resulting profile will be 
the convolution of the two Gaussian line profiles, i.e 



t(R) = Vf^therm(-R)2 + (Ttn^hiR)^ 



(E.7) 



In a planet-free disk, in which the gas parcels are moving 
on circularly Keplerian orbits, the line center at vo fundamental 
frequency shifts due to the Doppler shift is 



Av{R, 4>, i) = 



vo j GM, 



cos(0)sin(/). 



(E.8) 



Here we neglect that the massive planet and the host star are 
orbiting the common center of mass, instead the center of mass 
is set to the center of host star Moreover, the influence of the gas 
pressure on the angular velocity of gas parcels, which causes 
slightly sub-Keplerian orbital velocities due to radial pressure 
support, is also not taken into account. 
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